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ABSTRACT 


Но рБеһаулог Of the phase error probability 
density function of a first-order phase-locked loop in the 
presence of additive white Gaussian noise is determined 
using numerical techniques. In addition, the modulo-27 
probability density funetion and the cycle slipping density 
function are found. The algorithm used in these analyses 
involves numerical integration of a normalized and factored 
Fokker-Planck equation. Results are shown for cases in- 


volving various signal-to-noise ratios and initial conditions. 
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ДІ. INTRODUCTION 


A. BACKGROUND INFORMATION 

The space age spurred great interest in the development 
of phase-locking techniques since fixed tuned receivers were 
impractical for receiving signals experiencing a wide range 
οἱ doppler shifts. Due to the inherent nonlinearity of the 
loop phase comparator, early attempts to analyze phase- 
ESchodedooosbehavsropr involved the use of graphical phase 
plane methods which are summarized by Viterbi [1]. The 
initial use of Fokker-Planck techniques to determine the 
το "6 probability distribution of the first order 
loop phase error was accomplished by Tikhonov [2],[3]. Since 
then much work has been done using the Fokker-Planck method 
of analysis but little has been published dealing with the 
transient behavior of the phase-error process. 

Various techniques have been employed in attempting to 
Statistically describe the phase-error process of first 
order phase-locked loops. Stationary phase-error distribu- 
tions are well documented by Viterbi [1] and Tikhonov [2] 
and [3]. Transient solutions of the Fokker-Planck equation 
have been obtained by using eigenfunction expansions [4] 
and by exploiting the similarities between the Fokker- 
Planck equation and a one-dimensional heat-flow equation 
[5]. La Frieda [4] arrived at a statistical description 


of the reduced modulo-2n phase-error process by applying 





separation of variables. This application yielded an eigen- 
value problem which he solved numerically. Dominiak and 
Pickholtz [5] arrived at a direct numerical solution to 

the transient behavior of the phase-locked loop (PLL) by 
subjecting the Fokker-Planck equation to the same numerical 
procedures as was provided for the one-dimensional heat- 
flow equation by von Neumann and Richtmyer [6]. Their 
resulting partial difference equation was solved using a 
digital computer. 

This study was started in the belief that the previous 
Ori in the analysis of the transient behavior of the phase- 
locked loop could be expanded and improved. La Frieda's 
technique is applicable only for the reduced modulo-27 
solution and does not describe the actual behavior of the 
phase-error process. Dominiak and Pickholtz create serious 
doubt as to the completeness and validity of their technique 
ШОШ ОШОО снеге econleusion ОГ insignificant buildup of proba- 
bility mass at o = tk2T conflicts greatly with qualitative 
estimates of Viterbi [1] and Lindsey [7] and with one's 
Ion. ln was also felt that Dominiak and Pickholtz 
did not make a good choice of SNR for their numerical exam- 
peste sveartionary probability density function for their 
examples is so flat that any buildup of probability mass 
would probably be very small. 

In this paper the Fokker-Planck equation is normalized 
using several constants defined for the purpose [1] and [5]. 


The resulting equation is then expressed in terms of the 





product of two functions, one of which is indepen- 

dent of time. ПИК ОО Ше ЕГИП ООШ Т which is much 

шоо ега оле original probability density, is found 
numerically on a digital computer as a function of time. 
Merteıplying the two yields the solution to the Fokker- 
Planck equation which quantitatively describes the transient 
behavior of the phase-error process of a first order phase- 
locked loop. It is found that this procedure is very effi- 
cient computationally and readily provides numerical results 
describing the phase error process, the modulo-27 phase 


error process, and the cycle slipping statistics. 


B. REVIEW OF PHASE-LOCKED LOOPS 

A phase-locked loop is a device which automatically 
СОПОТ е an oscillator or periodic function generator so as 
to operate at a constant phase angle relative to a refer- 
ence signal source. Figure 1 shows the fundamental hard- 
ware structure of the loop with the appropriate signals. 
Bae three major components of the PLL are: 

a) A multiplier which multiplies the input signal to 
the loop and the output signal of the voltage controlled 
ео |аьор, 1.е., | 


кор) PAK.sins(t)cose'(t) 


Е 


AK, tsinle(t)-e'(t)] + sin[6(t) + 6'(t)]), 
| (1) 





where 0(t) and 8'(t) are the input signal phase and the 
voltage controlled oscillator output phase respectively. 
b) A linear low-pass filter with zero-input response 
eg (t) and impulse response f(t) which filters out the sum- 
frequency component of x(t) while passing a filtered version 
of the difference-frequency component, i.e., 
L 


e(t) = e (t) + х(%-ҒОдаА %>0. (2) 
0 


s ОШ сопсгоІ1ва о5с11Лабог (УСО), the output 
mequeneyyot which is equal to the sum of its quiescent 
ΙΙ, ТЕГЕН the product of its input voltage and a 
Beportienality constant, i.e., 

"ΤΙ Е 
БЕ = © + K,e(t). (3) 
Following through on the mathematics [1] it can be 


shown that for eS (t) - 0, 


£ 
f £(t-1)AK,sinle(1)-0' (1) Idd. 


5 (4) 


mE 
τ, EU. 


А 
For phase error ¢o(t) = 0(t)-0'(t) and total loop gain 


K = K. K. 


Ë 
EN ES f(t-A)sine(A)da3. (5) 
0 











Letting 





6, (t) - e(t) - в (6) 
0. (t) = 8g'(t) - o oU, (7) 
tnen 
9(t) = 0,(t) - 0,(t) (8) 
and 
de(t) de. (t) t 
Е = - AK Р J (S) Oo Ауа (9) 


which implies the PLL model shown in Fig. 2. 
Additive noise at the input which is narrowband Gaussian 
noise n(t) with zero mean and spectral density No near the 


signal frequency, can be expressed as 


all) -4J2tn, (t£)sinu,t + n,(t)cosuw_t) (10) 


where n, (t) and n. (t) are independent, zero mean Gaussian 

P ocesses of identical spectral densities, which are the 
same as the spectral density of n(t) but translated downward 
in frequency so as to be centered about zero frequency. A 
derivatıon similar to the derivation of (9) is made [1] with 


the noise having the form 


n'(t) = -n, (t)sins,(t) + n,(t)cose, (t) a) 


where n'(t) may be treated as a white Gaussian random 


process with spectral density Bor This leads to the model 
2 


10 





of the PLL in Figure 3, [1]. The equation of operation of 


the phase-locked loop is now given by 


de. (t) t 
_— Kf f(t-A)[Asiné(X) + n'(A)]dáàÀ. 
0 


τ 


EE) 
dt dt 


1i 
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FIGURE 1. Fundamental Hardware Structure of the 
Phase-Locked Loop. 


e(t) * FILTER 


de„(t) 
dt 
FIGURE 2. Phase-Locked Loop Modeled Without Noise. 
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FIGURE 3. Phase-Locked Loop Modeled With Additive Noise. 
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II.  FOKKER-PLANCK EQUATION 


Let a first-order phase-locked loop have an input signal 
e,(t) of constant frequency w,i.e., (t)=wt, and a voltage controlled 
oscillator with quiescent frequency ως" Additive thermal 
noise, n(t), at the front end of a receiver is a zero-mean 
Gaussian process whose spectral density is essentially flat 
over the frequency range of the receiver and may be assumed 
to be white with value "o. The equation describing the 


2 
Hoop Operation [1] is from (12) 





d$(t) = (u-w,) - KlAsing(t) + n'(t)] (13) 


where $(t) is the difference between the incoming signal 
phase and the VCO output phase and K is the total loop 
gain. 

It ean be shown[1] that e(t) is a continuous Markov process. 
The noise n'(t) is a stationary Gaussian process with zero 
mean and it is nearly white over the frequency range of 


interest. Therefore, 


N 
! ! = -9 - 
Е{п'(&,)п'(%)] = -” 8&(%„-®,). as) 
Under these conditions, the probability density function 
p(¢,t) of the Markov process phase error $(t) must satisfy 


the Fokker-Planek equation [1], i.e., 


13 





ὠμός. ὃ. [(w-w -AKsing )p(o,t)] 


ət 9ф 
KÊN 2 
O οφ 5) 
SIC 5 


9ф 
mtc initial condition 


р(%,0) = 6(Ф-Ф) 


and normalizing condition 


со 


/ р(%,5024% - 1 


— © 


It is convenient to express (16) in normalized form. 


В, = AK/4 loop-noise bandwidth 
« = A*/N B, signal to noise ratio 
у = (0-0 )/4B, frequency ratio 
т = ШВ,% normalized time. 


L 


Substituting (19) - (22) into (16) yields a normalized 


Fokker-Planck equation which is 


2 
Ho E ag [(siné-y)p(6,1)] t 1 3 pesar) 
дф 


An analytic solution to (23) has never been found. 


(16) 


CIT) 


(18) 


Let 


(19) 
(20) 


(21) 


{66} 


оз 


It is apparent that by considering only very small phase 


errors, l.e., sine = Ф, then (13) can be linearized. 


14 


Making 





the same approximation in (23) gives a closed form solution 
which satisfies the initial and boundary eonditions on p(6,:). 


fiers solution is 


я L a[#-v-($ -v)e J* 
р(ф,т) = м ee a Е (24) 
ие <‘) 2(1-е * ) 

which is Gaussian with 

mean = m(t) = y + (Ф, ΕΣ: (25) 
and 

: _ 2 _ =21 
variance = с (т) = (1-е )/a . (26) 


It is the Gaussian nature of (24) which provides the 
information necessary to choose the initial boundary condi- 
Benzsrunetion Гог the numerical algorithm used in the solu- 
tion of (23). However, before the selection of an initial 
ιν ls discussed, it is appropriate to analyze the 
етпе моо сеа in arriving at the numerical algorithm 


used in this study. 


1 





III. TRANSIENT ANALYSIS 


А. FUNDAMENTAL TECHNIQUE 
The fundamental technique employed in this study is the 
treatment of the normalized solution to the Fokker-Planck 


коор equation as the product of two functions, i.e., 


ВЕ = 2(9,:) ехрГоа(совф + ү9)1. (27) 


This exponential factor was selected on both mathematical 
and intuitive grounds.  Intuitively, one would expect that 
after having slipped a cycle the loop dynamies would again 
exhibit some measure of steady state behavior prior to 
another slip. Lindsey [7] makes a similar observation 
in his analysis of cycle-slipping probabilities. At any 
one time the loop phase-error tends to wobble about one of 
an infinite set of points spaced by 2v radians. It seems 
reasonable that the phase error probability density should 
resemble a mod-2n stationary density function in each en 
interval with a broad taper away from the origin reflecting 
τ АТК i.e. cycle slipping, nature of the problem. 
The "envelope" encompassing this succession of quasi- 
mod-2n functions can be easily imagined to be much broader 
and smoother than the p(¢,t) function. It seems plausible 
that the function g(¢,1) from (27) is this smooth "envelope". 
The desirability of this observation is that the "smoothness" 


ОЛИ ООП would allow a much coarser grid to be used in a 





numerical integration scheme than would be possible with 
the function р(%,т) 
Applying (27) to the normalized F-P equation (23) one 


obtains 


2 
2609.1) = (ү-віпф) — + 12 209.1) 2 (28) 
9% 


It is now apparent that ООШ ШЕ cor eto.) is eomputationally 
simpler than solving for p(é6,:). And, since the exponential 
factor is a constant for each value of » along the $-axis, 
р(ф%.т) can be reconstructed readily using (27). 

Equation (28) is evaluated numerically by using the 
central-difference method to approximate the derivatives, 
such that with Ad and At being the increments in $ and τ 


respectively, then 


2600.1) 2 glesrtaricele,r) (29) 
3g(6,1) . g(ot^o,1)-g(0,1) (30) 
дф Аф 


3΄ρ(ὀ;τ) .γΕ(ὀ1Δὀ,τ)-β(ὀ»τ) βίϑ,τ)-βίὀ-λὀ;τ)η 1. 
9 Аф 8% ab 


дф 
.Е(Ф%Аф,т) - 26(%,т) % р(9-Аф,т) 
x _ (ST) 
(4%) 


and it follows directly that 


(32) 
509 ,т+Ат)=609,т)+(у-віпф) Cg (646 :)-g( 6746, 1)] 


τ 


π---πίείφεδὀ.,τ)-2β(ϕ.,τ)ερίϕ-δϕ»τ)]. 
) 


1 





B. OTHER CONSIDERATIONS 

iL) S] Condition 

It is the response of (23) to an input impulse 

function (17) which is sought. However, it was felt that 
а finite difference method algorithm acting on the extreme 
ο. ο ΙΙ ντες θα numerical delta function of Ad width and 
1/49 height might result in significant errors in the re- 
sults. Therefore the Gaussian nature of (24) was utilized 
επ initial boundary condition for the algorithm. 
The boundary condition selected was a Gaussian function 
which corresponds very nearly to the response of (23) at 
some very small time, assuming the impulse function (17) 
existed at time т = O. To check the accuracy of this boun- 
dary condition three techniques were employed. The first 
was to plot the variance of the calculated phase error 
versus time to determine whether or not the loop is opera- 
Ling in its linear region; if so, the variance vs. time 
curve would extrapolate linearly back to an equivalent 
origin, where the variance is zero at т = 0. The second 
was to actually test the phase error response to the numer- 
1са1 delta function mentioned above. The third was to 
compare stationary and transient mod-2- results with those 
suc3dpubdédshed [3], [4] and [5]. 

When the initial condition is the Gaussian function, the 
equivalent time, το ι.. επεπις.τοπεαρ”οη would, occur 


if p(o,0) = $(670,) is given by (26) as 


MS 





cog UG (33) 


and it follows that 


JE Г1п(1-ас 277. (31) 


chars study B: t = 0 anda B: t = E are the times at which 
the equivalent delta function and the Gaussian initial 
condition existed respectively. In [5], however, B; t = 0 
was taken to be the time of existence of their Gaussian 
um πο eonda cion, and this does not truly represent the 
dependency of p(6,1) on time. 
2. End-value Algorithm 

The use of a digital computer in evaluating this 
numerical integration scheme requires a finite interval 
which comprises the sample points along the ф-ахіѕ. А 
ШЕ τάσδε source of error in any computations per- 
formed on this fixed interval is the method by which the 
end points of the interval are evaluated after each At step. 
The central-difference method algorithm in this study, if 


applied exclusively, would result in an error propagating 


ο ο οπς ие гешепс рег a, increment starting with the 


first лт increment. This error would make the procedure 
useless for large t unless the total ¿-interval were made 


prohibitively large. To minimize error from this source, 


19 





a sixth-order Stirling extrapolation formula is used to 
evaluate both of the end points. | 
ο "ουσ. τον 

Stability in a numerical algorithm implies that the 
errors generated at each step of the numerical integration 
process tend to decay as the solution progresses in time. 
The method used in this study to determine the stability 
conditions is patterned directly after von Neumann's Fourier 
series technique [6] using the format presented in [8]. The 
mesh is assumed sufficiently small such that the coeffi- 
cients in (32) may be treated as constants over a region 
which is small but contains many mesh points. For N sample 
points spaced Ad radians along the ¢-interval, the initial 
errors, i.e. t = 0, which are assumed to be located at the 
Sample points can be expressed for each point in complex 


exponential form as 


it 28 К(4%) 
r ο 


E. = у An Ра ο ο ο Ne] (35) 
n=0 
where 
mon 
а = ЕО C ο 


Because the finite-difference equations will always be 
linear and, therefore, separate solutions additive, one need 
only consider the propagation of the error due to a single 


е18К(А%). 


кект E.g.: The coefficient A is a constant and 


20 





соп ше neglected; and for simplicity the subscript п оп 
B is dropped. 
enve Sti gate the propagation of this error for 


increasing t, one must find a solution of the finite-diff- 


erence equation which reduces to οὐ ΒΚ{ΔΦ) when t = &(At) 
Tg is defined as the error at the us point along the ф 
em 


шәрегуат ас опе L£ iteration along т. Assume that the 


error has the form 


εὐ)θκίδϕ) ¿AT . ,J8k(A0) ¿AL(A7) 
k2 


ET: (37) 


DP ο ο Бепсогас; і5 a complex constant. 


јвк(лф) 


uere £ = е 
Clearly, (37) reduces to e when £ = 0. Also, the 


1] τοῖ increase Tor increasing r provided 


τ. (38) 
Since the error Pro satisfies the same difference 
tons taste, 7), substitution of (37) into (32) and 
j8k(^9) ,2 


subsequent division by e A e tas 


m 
! 


1 +(y- sind) ро exp[-j8(49) Ј?} 





+ і 
а 


(exp[j8(^9)1-2 * ехр[-јв(лф) 1} 
А 


1+2 (y- sin9)Htsinle (ag) TÉ y —5 (cos[g(a9))- -1}. 
(29) (39) 


za 








If one assumes that the 2(y-sine)gasin[8(a9)] term is 
negligible for large N and very small. Ad, then from (38) 


the following inequality is established: 
a 2 
Ат < 5(Δ9) : (40) 
For numerical calculations let At and Aq be related as 
a 2 
Ат = Τ(Δ9) 2 (H1) 


ОО Лес ОГ о in the vicinity of unity and for Ap << 1, 
рлеп Іт << А9. Hence, the dropping of the term containing 
E унес σα. 
4. Cycle Slipping 

En Sono tab sorbine boundaries around some 
sub-interval of $, e.g. [--,-7], enables one to determine the 
probability that the phase error has not exceeded the limits 
ο. οπου sub-interval at any point in time. As the density 
ото spreads out, any probability mass which reaches 
the limits of the sub-interval is removed from the system 
ештегеіу by forcing the density function кере жене outside 
Teemu interval. The effect of this technique is to make 
the remaining area of the density function equal to the 
probability that the phase error has not reached the limits 
of that sub-interval. 


The cycle-slipping density function is denoted 


ҮЙЕ (55 distinguish it from the corresponding ОЕ ол 





р(ф,т) for the unbounded case. As long as the phase error of 
ШЕР 1σορ remains within the limits of a Specs Sula етта. 


the density function q($,1), with initial condition 


q($,0) = 6(¢), will satisfy the Fokker-Planck equation (23), 


т.е. 
2 
aqlo,r)_ ale) + l3 q(0,1) — (uo) 
дт 9 $ α ТБ | 
Following the development from Viterbi [1], for y = 4 = O 


© 
the probability that the phase error $ has not reached some 


value Ф, O MeT is given Dy 


Ww 
y(t) = εἴ а(%,т)94ф (83) 
2, 


where 


q(6,, 1) = q(-$,,1) = 0 for all t . (44) 


The fundamental difference between q(6,1) and p (6,1) is 
МИРЕ“ 409%9,т) іс пос strictly a probability density function, 
unless it is normalized by w(t). The numerical algorithm 
used to generate q(¢,1) is the same as for generating р(%,т) 
except that the end points ale T) and ο το ут) are forced 
to be zero. The end points for p(6,1) are found by extrapo- 
Шом сапе simulating an infinite interval. 

In this study an analysis was made of the cycle 


A roba Dilities ror the primary cycle slip, i.e., 


[-¢, 59, ] = [-т,т] (45) 


ББ 





anderorscher primary half-cycle slip, l.e., 
[-¢,59,] ж (-т/2,т/2| ° (H6) 


Plots of ¥(t) versus t show the time dependent behavior of 


ems cycle slipping probability. 
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IV RULES 


The entire numerical analysis of this study was pro- 
grammed for an IBM-360/67 computer using Fortran IV language, 
level G. Doubie precision computation was used exclusively 
Са шаге round—off and truncation errors. The initial 
condition on p($,1) in (23) was a Gaussian function with 
mean $, and standard deviation a "7/20. The $-interval of 
investigation was [-11-,11-] with g($,1) and, therefore, 
p(6,1) being computed over the entire interval to allow for 
the non-symmetrical conditions Y Z 0 and e о: 
considerations specify that Ar < ar and to accom- 
mođate this requirement Ат was chosen such that 
Ат = ο... The sampling interval along the d-axis 
was selected to be Ad = 1/50 except in the cycle slipping 
calculation where Ad = 7/100. For the examples which will 
be shown here, it was found that negligible improvements 
were obtained for smaller values of ^4 and со: ue acce 
a "delta function" of width 4$ and height 1/49 was used as 
the initial condition on p($,1) in several cases, and the 
results obtained varied only by approximately 1 part in 
20,000 from those obtained with the initial Gaussian condi- 
tion above. Table 1 indicates the various conditions which 
were studied. The specific values were chosen so as to 
provide comparison with previously published work [1], [4], 


(51. Also included in the table are the accuracy 


29 





considerations for this study expressed as areas of the 
Шыр density function aiter the indicated number of 


LERNEN Steps, 


TABLE T 


Number of Area Under 
Ат Time Steps Curves 
ER = 20,000 1.000 
ws 10,000 1.000 
5η, 1) 20.000 1.000 










Figures 4, 5, 6, and 7 are graphs relating to the condi- 


pcs = 1-07 ап4 ү - = 0. Figure H is a sample of the 
dynamics of the g(#,1) function from (27). It demonstrates 
the "smoothness" of the function with respect to the р(ф,т) 
Mco ch is shown in Figure 5. Transitional and 
stationary modulo-2- densities are depicted in Figure 6. 

IN cur. mosulo-es density function is seen to be a 

e odneo of the curve determined by Viterbi [1]. Of 
particular interest from Figures 5 and 6 is that the р(%,т) 
ιτ continues to spread out with buildups of probability 
mass at the tk2n points even though the steady state modulo-2T 
density function has been reached. Figure Y is a plot of the 
variance of the phase error probability density function 
versus time. It is interesting to note the linear nature of 
επ ο νε τος Br t 6. 
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La Frieda's eigenfunction analysis [4] resulted in 
transitional and steady state mod-27 densities which were 
reproduced using the algorithm of this study. Figure 8 
is included to show the p(o,t) curves corresponding to 
each of the times used in [4] with an additional curve 
to demonstrate clearly the asymmetric probability mass 
buildup associated with the conditions a = 2, y = 0, апа 
$6 = T/2. 

Inasmuch as La Frieda used a different normalizing 
scheme, Dr, in his study, equivalent times, i.e. B t=Dr/2, 
are shown on 'the graph for ease of reference. The 
asymmetry of probability mass distribution is due to the 
main body of p(¢,1) being offset from $ = 0 for some finite 
КООШ ОГ Ое while It is shifting. Since р(ф,т) 15 
settling even while it is shifting, more mass will accumulate 
at ¢ = 2n than at y = -27. Once the main body of p(¢,1) 15 
settled about ¢ = 0, however, the slippage of mass in 
either direction becomes equally likely and the level of 
j mm ry |! remain constant for increasing time. 

Viterbi dealt with a detuned loop which was also ana- 
lyzed in this study. For this case a = 1.0, y = sin(1/4), 
and 2 8 О [1]. Figures 9 and 10 show the mod-21 prob- 
ability density function and p(94,71) for these conditions 
respectively. As before the steady state mod-2r density 


funetion determined by this study is a reproduction of that 


found by Viterbi. As would be expected and is demonstrated 
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ШШШ ШО спгуез the detuned loop exhibits a definite 
tendency to slip in the direction of the steady-state phase 
shift caused by detuning. 

ПЕЕ is ἔπε ΡίΦ.τ) funetüon corresponding to 
« = 0.1 and y = $ ου πας бет ОГ conditions is one of 
those selected by Dominiak and Pickholtz [5] in their 
analysis of the first order PLL by a numerical technique. 
This curve shows the buildup of probability mass which 
went undetected in the other study. It is felt that Domin- 
jak and Pickholtz did not find this buildup for several 
reasons. Firstly, they did not look look enough in their 
time frame. Secondly, their attempt to deal directly with 
the p(6,1) function rather than an intermediate function, 
e.g. g(¢,1) in this study, resulted in a grid mesh size 
which, although mathematically correct, was too coarse to 
detect the subtle and yet significant mass slippage of the 
Paes-ipaelıvy density function. Finally, From rl ats 
СОР ЛСП ШЕ АГ tor a = 0.1 the mod-2r7 density function is 
МОО О ОШОО пс mot a good choice of SNR with which to 
ου οι ρποσαυα τιν mass buildup. 

Figures 12 and 13 show the dynamics of the cycle slip- 
Papa dencesarunetions, d(¢,1), for primary cycle slips, 
[-7,7], and primary half-cycle slips, [-1/2,7/2]. These 
τ ος Ροτεα Гог the conditions a = 0.5 and $o = ү 50. 
In Figures 14 and 15 are the curves which demonstrate how 


οι j Tipping probability, v(1), varies with time for 
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the cases shown in Figures 13 and 14. The function а(%,т) 
ШОШ ОНЫ шлу оо е slip is essentially identical to the 
Ho ООЛО ОП ксі! ә significant amount of probability 
mass is lost at the absorbing boundaries at 4 = іп. Until 
this loss occurs the area under the q(¢,t) curve, i.e. v(t), 
is 1.0. The time at which y(t) drops below some value, 

e.g. 0.5, can be determined from Figure 14. From Figures 

13 and 15 similar data may be obtained for the primary 


half-cycle slip which has absorbing boundaries at ¢ = + = 5 
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V. CONCLUSIONS 


This study has resulted in the development of a valid 
technique for numerically finding transient solutions for 
the Fokker-Planck equation as applied to the first-order 
PLL. That the technique works is demonstrated by its repro- 
duetion of existing published stationary mod-2- solutions 
in [1] and transient mod-27 solutions in [4]. Results 
not previously available for the phase-error probability 
density function and the cycle Slipping probability are 
presented here. Most important, however, is that the 
ο ο ЛЕУ е Notion or the buildup of probability mass 


around points spaced 27 radians has been verified. 
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